Supplemental Information 

A Relationships between enzyme masses and flux control coefficients 

According to Brown [12], in the state of minimum total enzyme concentrations, the flux control coefficient 
of any enzyme will be proportional to that enzyme concentration for linear metabolic pathways. This 
relationship is given in Equation (SI) where Cf represents a scaled flux control coefficient. 

1 = c (SI) 

To give more insight in this relationship, the derivation of this equation will be explained with the two-enzyme 
metabolic network shown in Figure SI. 
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Figure SI: Simplest biochemical pathway with two enzymes {E\ and E 2 ), two fixed metabolites (S and 
P), and one floating metabolite (X). Floating and fixed species are indicated in blue and red squares, 
respectively. Fixed species are underlined. 

Reducing E\ by an infinitesimal amount (5Ei) results in a pathway flux decrease. An increase of E 2 by 
an infinitesimal amount (8E 2 ) is necessary to return the flux to its original rate. For this model, we assume 
the model is in the state of minimum total enzyme concentrations for a particular flux. In other words, the 
minimum total enzyme concentration is found for a particular flux, which means that Equation (S2) must 
hold. 

8 Ex + 8E 2 = 0. (S2) 

Otherwise, enzyme concentrations can be adjusted such that the total enzyme concentrations decreases and 
the flux remains constant. Because 8 J — 0, the relationship in Equation (S3) must be true. 

SJ C( -8Ex + ^-8E 2 = 0. (S3) 



J Ex E; 
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Rearranging Equation (S3), in observation of Equation (S2) gives Equation (S4). This Equation can be 
generalized into the earlier given Equation (SI). 



hiX Ej-x 

However, Klipp and Heinrich [13] show that this relationship is not valid for branched pathways, and derived 
a new relationship between flux control coefficients and enzyme concentrations in the state of minimum total 
enzyme concentration. This relationship is given in Equation (S5). Here, C and e represent the flux control 
coefficients matrix and the enzyme concentration vector, respectively. See [13] for a proof of this relationship. 

C T e = e. (S5) 
From this equation, Equation (SI) can be derived as follows: 



C'x 1 C^\ T (Ex\_(Ex 



C( 2 Ci 2 J \E 2 J \E; 

C J l 1 Ex+C J l *E 2 = Ex 
C^ 2 Ex + C'2' ' E 2 — E 2 . 



(S6) 
(S7) 
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Since Ji = J 2 in the linear pathways, the above equations can be simplified: 

Cf(E 1 +E 2 ) = E 1 



Ci(E 1+ E 2 )=E 2 (S8) 



(J J _ ij 



J| (S9) 
where 

^T = ^^i- (S10) 

Metabolic cost of enzymes can, however, vary due to their different molar masses. As a consequence, 
weight factors need to be introduced into Equation (S10). This relationship is given in Equation (Sll), 
where m; is the molar mass (g x mole^ 1 ) of E i7 and Mt the total enzyme mass. 

M T = '£ i E i m i = Y,M i . ( sn ) 

The addition of the molar masses to Equation (S10) results in the modification of Equations (SI) and (S5). 
We illustrate this with the linear model shown in Figure SI. Reducing E\ by an infinitesimal amount (5Ei) 
results in an increase of E 2 by an infinitesimal amount (SE 2 ). In an optimized state — minimum total 
enzyme mass for a particular flux — Equation (S12) must hold: 

6M T = SE imi + 5E 2 m 2 = 0. (S12) 

Moreover, SJ = 0 in the optimized case, which is given by Equation (S3). Rearranging Equation (S3), in 
observation of Equation (S12) gives Equation (S13). 

Cl - ° 2 (S13) 



Eim\ E 2 m 2 

Next, Equation (S13) can be generalized to Equation (S14). This relationship shows that, in the state of 
minimum total enzyme concentrations for a given flux, the flux control coefficient of any enzyme will be 
proportional to that enzyme mass. 

p— - c. (S14) 

This relationship is similar to Equation (SI), whereas it incorporates the effect of molar mass. A similar 
adjustment was done for Equation (S5): 

C T e m = em, (S15) 

where the z-th element of em is Eiirii. 

From Equation (S15), one can derive Equation (S14) for linear pathways as follows: 



(C^ Cf\ T (E 1 \(m 1 \_{E 1 



cl-) &){Z)-\$)Q) <S16) 

C^E imi + C( 2 E 2 m 2 = E imi ( 
C J 2 2 E imi + Ci 2 E 2 m 2 = E 2 m 2 . [ ' 



Since Ji = J 2 in linear pathways, 



C((Eimi + E 2 m 2 ) = Eimi 
C 2 (Eimi + E 2 m 2 ) = E 2 m 2 



ru] E\ra\ 

r<J E 2 m2 

°2 - Tf7 

cj = Cj = _j_ 



(S18) 

(S19) 
(S20) 



26 



B Maximize metabolic flux given a total enzyme mass constraint 

In this section, we propose an optimization algorithm, which was used to optimize steady-state flux given a 
finite total enzyme mass. In such an optimized state, according to Equation (S14), flux control coefficients will 
be proportional to enzyme mass. We developed an algorithm that finds, independent from initial conditions, 
the highest flux through a metabolic pathway. This algorithm iteratively increases enzyme concentrations 
with low C/ fe x (Einii)^ 1 and decreases enzyme concentrations with high C/ fc x (.EjTOj) -1 . Here, k represents 
the flux to optimize. Eventually, C*/ fc x (Eiin^ 1 converges for each enzyme in the branch to optimize, such 
that Equation (S14) holds when this algorithm is finished redistributing enzyme concentrations. 

Note that enzyme concentrations can tend to zero if the network has, for instance, more than one 
independent flux. Because this algorithm makes percentage changes to enzyme concentrations, enzyme 
concentrations can go asymptotically to zero. For this reason, only C/ fc x (Eirrii)^ 1 are taken into account 
if the enzyme concentration is greater than a threshold close to zero. Enzymes favorable in optimized states 
can be close to zero in non-optimized states. Therefore, an exception is made for enzymes that start with 
an enzyme concentration less than the selected threshold, which they lose if they start to decrease. 

The first step of this algorithm is to perform steady-state analysis. Secondly, scaled flux control coef- 
ficients are determined and divided by the current enzyme concentrations times the corresponding molar 
mass. Subsequently, the algorithm proposes to increase the enzyme with the highest C/ fc x (Eiirii)^ 1 and 
decrease the enzyme with the lowest C/ fc x (Eiirii)^ 1 . Proposed changes in enzyme concentration are based 
on the enzyme with the lowest Eirrn. Moreover, these proposed changes must be sufficiently small, because 
flux control coefficients alter when enzyme concentrations are altered. Because of the constraint given in 
Equation (Sll), Equation (S21) must hold while redistributing enzyme concentrations. 

SM T = m max 5E max + m min SE min = 0 (S21) 

Here, E max and E m i n represent the enzymes with the maximum and minimum C- 7fc x {E.mi)- 1 where E t 
is greater than the threshold. Hence, the proposed change for the enzyme with the lowest C/ fc x {Eim^ 1 
is given by Equation (S22). 

^■^max ^min (S22) 

During every iteration the concentration of two enzymes changes, while the total enzyme mass remains 
constant. The system, subsequently, goes to a new steady-state with a corresponding flux. The optimal 
enzyme distribution is found once the flux to optimize starts to decrease. An example of this numerical 
algorithm is shown in Figure S2. 
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Figure S2: Simulation results for the maximization of pathway flux for a six-enzyme linear pathway. 
Increasing pathway flux, (b) Converging C/ fe x (.E^m*) -1 . 
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Figure S3: Branched metabolic network. Case 1: no feed-back, case 2: feed-back from X 2 to Vi, and case 3: 
feed-back from X4 to V\. 



C New relationships between flux control coefficients and enzyme concentrations for 
branched networks 

This section describes numerical results supporting new relationships between flux control coefficients and 
enzyme concentrations similar to the Brown's result, Equation (SI). Our result presented here is funda- 
mentally different from the Klipp's result in that it provides extra information than Equation (S5) can 
do. 

First, consider case 1 of the branch network described in Figure S3. We minimized total enzyme mass 
under the constraint of fixed fluxes. Figure S4 illustrates the system's behavior. Figure S4(a) confirms 
that the fluxes remained constant during the simulation. Figure S4(b) shows that some enzyme levels 
increased, while others decreased during the minimization of total enzyme mass. Figures S4(c)-(f) illustrate 
C^ h {EiTni)~ 1 values during the optimization procedure. As can be seen from these figures, C i k (Eimi)^ 1 
values converged for linear segments in branched pathways, just as they do in linear pathways as was shown 
previously by Brown [12]. These results suggest that Equation (S14) also holds for linear segments in 
branched pathways. 
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(S23) 



We performed simulations (not presented here) for larger branched networks up to twenty metabolites, 
where the model contained more than two branches. Results, again, indicated that Equation (S14) holds 
for linear segments in branched pathways in the state of minimal enzyme mass for a particular set of fluxes. 
This suggests that Equation (S14) can be modified into Equation (S24). 

(S24) 



Eivrii 



if k and i are in the same branch. 

We further investigated the network shown in Figure S3, but now with different types of feed-back (e.g. 
case II and III). These results are illustrated in Figures S5 and S6. Here, C/ fc (Eirrii)^ 1 did not convergence 
for branches where the negative feed-back originated, which suggests that Equation (S24) need to be adjusted 
into Equation (S25). 

-^- = c, (S25) 
if k and i are in the same branch without an originating feed-back. 
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Figure S4: Simulation results for the minimization of total enzyme concentration while fixing fluxes for a 
six-enzyme branched pathway without any feed-back (Figure S3. The C/ fc (Eiin^ 1 of reaction sets from 
different branches — v± and v<2, V3 and v§, and and vq — did converge. For parameters, see Section F.4.2. 



30 





0.8 




0.6 








0.4 


o 


0.2 




0 




-i = 1 
■ i = 2 

i = 3 
■i = 4 
■i = 5 

i = 6 





0.6 




0.4 


1 


0.2 


6 






0 


O 


-0.2 




-0.4 




-0.6 



0 500 1,000 1,500 2,000 
Iterations 
(a) Convergence of C/ 1 (Eiim)^ 1 



0 500 1,000 1,500 2,000 
Iterations 
(b) Convergence of C/ 2 (Emu)' 1 




500 1,000 1,500 2,000 
Iterations 
(c) Convergence of C/ 3 (Eiirii)^ 1 



500 1,000 1,500 2,000 
Iterations 
(d) Convergence 




0 500 1,000 1,500 2,000 

Iterations 
(e) Convergence of C/ 5 (Eirrn)^ 1 





0.2 - 
0.1 - 
0 



■ *tv'v. >i' V' « , , 't>V' ' l ' | h 



= 1 
= 2 
= 3 
= 4 
= 5 
■-- i = 6 



0 500 1,000 1,500 2,000 
Iterations 
(f) Convergence of C(* (Enm)' 1 



Figure S5: Simulation results for the minimization of total enzyme concentration while fixing fluxes for 
a six-enzyme branched pathway with negative feed-back from X 2 to v\ (Figure S3. The Cf h {E i m i )~ 1 of 
reaction sets from different branches — V\ and v 2l v 3 and Ug, and V4 and Vq — did converge. For parameters, 
see Section F.4.2. 
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Figure S6: Simulation results for the minimization of total enzyme concentration algorithm for a six-enzyme 
branched pathway with negative feed-back from X4 to v\ (Figure S3). The C i k (Eiin^ 1 of reaction set V4 
and vq did not converge. For parameters, see Section F.4.2. 
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D Distribution of the deviation from the desired input-output characteristics 

The fitness of each GRN was based on the deviation from the desired input-output characteristics. Specifi- 
cally, the deviation was quantified by the relative difference between the desired and training output values. 
The distribution of these relative deviation values within the trained region was significantly different from 
a normal distribution (Shapiro- Wilk test [?] showed p <C 1.2 x 1CP 5 ). Quantile-Quantile plots given in Fig- 
ure S7 show that the distributions are rather similar to a gamma distribution, thus the standard deviations 
were determined based on a gamma distribution. 
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Figure S7: Quantile-Quantile (Q-Q) plots of the relative deviation values within the trained region against 
gamma and normal distributions. Shape parameter k equals 2.25, 2.25, 5.75, and 7.0 respectively. 
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E Perturbing concentrations of boundary metabolites 

In this section, we illustrate how large the differences between optimal enzyme concentrations can be for 
different concentrations of boundary metabolites. During training of GRNs the concentrations of boundary 
metabolite S and P were often between of 1 - 100 \xM and 0.01 - 1 fj,M, respectively. In Figure S8 we 
show the optimal enzyme concentration given a total enzyme mass constraint for two particular enzymes. 
Figure S8(a) shows that perturbation of S has a large effect on the optimal enzyme concentration, while 
perturbing P has almost no effect. In contrast, Figure S8(b) shows that both S and P have a significant 
effect on the optimal enzyme concentration of this particular enzyme. 
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Figure S8: Relationships between concentrations of boundary metabolites S and P and two particular 
enzymes, (a) Perturbing S has a much larger effect than P on the optimal enzyme concentration of this 
enzyme, (b) Perturbing both S and P has a large effect on the optimal enzyme concentration of this 
enzyme. 
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F Models of metabolic and gene regulatory networks 

F.l Two-enzyme metabolic network 

# PySCeS input file 

# Easy Pathway 

FIX: S P 
Rl: 

S = X 

El*kcatl*(S - X/Keql)/(S + KmS*(l + X/KmX)) 

R2: 

X = P 

E2*kcat2*(X - P/Keq2)/(X + KmX*(l + P/KmP)) 



# InitExt 
S = 75 

P = 25 

# InitPar 
El = 1 
E2 = 1 
ml = 1 
m2 = 1 

kcatl = 1000 
kcat2 = 500 

Keql = 5 
Keq2 =22.5 

KmS = 1 
KmX = 3 
KmP = 2 

# InitVar 
X = 1 
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F.2 Three-enzyme metabolic network 

FIX: S P 

Rl: 

S = X 

El*kcatl*(S - X/Keql)/(S + KmS*(l + X/KmX)+ alpha*Y/KmY) 

R2: 

X = Y 

E2*kcat2*(X - Y/Keq2)/(X + KmX*(l + Y/KmY) ) 

R3: 

Y = p 

E3*kcat3*(Y - P/Keq3)/(Y + KmY*(l + P/KmP)) 

# InitExt 
S = 10 

P = 0.01 

# InitPar 
El = 1 
E2 = 1 
E3 = 1 
ml = 1 

m2 = 1 
m3 = 1 

kcatl = 500 
kcat2 = 50 
kcat3 = 100 
alpha = 3 

KmS = 1 
KmX = 0.5 
KmY = 1.5 
KmP = 1 
Keql = 20 
Keq2 = 15 
Keq3 = 50 

# InitVar 
X = 1 

Y = 1 
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F.3 Five-enzyme metabolic network 

FIX: S P 



Xl/Keql)/(S + KmS*(l + Xl/KmX)) 



X2/Keq2)/(X1 + KmX*(l + X2/KmX)) 



X3/Keq3)/(X2 + KmX*(l + X3/KmX) ) 



X4/Keq4)/(X3 + KmX*(l + X4/KmX) ) 



P/Keq5)/(X4 + KmX*(l + P/KmP)) 

# Fixed Species 
S = 1 

P = 1 

# Variable Species 
XI = 1 

X2 = 1 
X3 = 1 
X4 = 1 

El = 0.7 
E2 = 1.3 
E3 = 1 
E4 = 1 
E5 = 1 



# Parameters 


kcatl 


= 100 


kcat2 


= 50 


kcat3 


= 100 


kcat4 


= 75 


kcat5 


= 125 


KmS = 


5 


KmX = 


0.5 


KmP = 


10 


Keql = 


-■ 20 


Keq2 = 


■■ 15 


Keq3 = 


= 10 


Keq4 = 


= 20 


Keq5 = 


= 15 



S = XI 
El*kcatl*(S - 

XI = X2 
E2*kcat2*(Xl- 

X2 = X3 

E3*kcat3*(X2 - 
X3 = X4 

E4*kcat4*(X3 - 
X4 = P 

E5*kcat5*(X4 - 
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F.4 Six-enzyme metabolic networks 

F.4.I Linear with feed-back 
FIX: S P 

Rl: 

S = XI 

El*kcatl*(S - Xl/Keql)/(S + KmS*(l + Xl/KmXl)) 

R2: 

XI = X2 

E2*kcat2*(Xl- X2/Keq2)/(X1 + KmXl*(l + X2/KmX2)) 

R3: 

X2 = X3 

( (E3*kcat3*X2/KmX2) * (1-X3/ (X2*Keq3) ) * (X2/KmX2 + X3/KmX3) ** (n-1) ) / 
(((X2/KmX2+ X3/KmX3)**n) + (l+(X5/KmX5)**n)/(l+alpha*(X5/KmX5)**n)) 

R4: 

X3 = X4 

E4*kcat4*(X3 - X4/Keq4)/(X3 + KmX3*(l + X4/KmX4)) 

R5: 

X4 = X5 

E5*kcat5*(X4 - X5/Keq5)/(X4 + KmX4*(l + X5/KmX5)) 

R6: 

X5 = P 

E6*kcat6*(X5 - P/Keq6)/(X5 + KmX5*(l + P/KmP)) 

# Fixed Species 
S = 100 

P = 0.5 

# Variable Species 
El = 1 

E2 = 1 
E3 = 1 
E4 = 1 
E5 = 1 

E6 = 1 
XI = 1 
X2 = 1 
X3 = 1 
X4 = 1 
X5 = 1 

# Parameters 
kcatl = 70 
kcat2 = 50 
kcat3 = 100 
kcat4 = 75 
kcat5 = 125 
kcat6 = 150 
KmS = 5 
KmXl =0.5 
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KmX2 =1.5 
KmX3 = 4 
KmX4 = 2 
KmX5 =0.75 
KmP = 7 
Keql = 20 
Keq2 = 15 
Keq3 = 10 
Keq4 = 20 
Keq5 = 15 
Keq6 = 30 
n = 1.5 
alpha = 1 



F.4.2 Branched 
FIX: S PI P2 

Rl: 

S = xl 

((El*klf*S/Ks)*(l-xl/(S*Keql))*(S/Ks + xl/Kxl))/ 

(((S/Ks + xl/Kxl) ) + (alphal*(x4/Kx5))+(alpha2*(x3/Kx5))) 

R2: 

xl = x2 

E2*k2f*(xl - x2/Keq2)/(xl + Kxl*(l + x2/Kx2)) 

R3: 

x2 = x3 

E3*k3f*(x2 - x3/Keq3)/(x2 + Kx2*(l + x3/Kx3)) 

R4: 

x2 = x4 

E4*k4f*(x2 - x4/Keq4)/(x2 + Kx2*(l + x4/Kx4)) 

R5: 

x3 = PI 

E5*k5f*(x3 - Pl/Keq5)/(x3 + Kx3*(l + Pl/Kpl)) 

R6: 

x4 = P2 

E6*k6f*(x4 - P2/Keq6)/(x4 + Kx4*(l + P2/Kp2)) 

# InitExt 
S = 10.0 
PI = 1.0 
P2 = 1.0 

# InitPar 
El = 1.0 
E2 = 1.0 
E3 = 1.0 
E4 = 1.0 
E5 = 1.0 
E6 = 1.0 
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klf = 1 
k2f = 1.5 
k3f = 1.0 
k4f = 1.0 
k5f = 2.6 
k6f = 1.0 



Keql 




0. 


7 


Keq2 




1. 


1 


Keq3 




2. 


1 


Keq4 




1. 


4 


Keq5 




0. 


4 


Keq6 




1. 


6 


ml = 


1 


,2 




m2 = 


0. 


.7 




m3 = 


2 


.7 




m4 = 


0 


,8 




m5 = 


1, 


.4 




m6 = 


1. 


,0 





Ks = 10 
Kxl = 5 
Kx2 = 2 
Kx3 = 0.5 
Kx4 = 4 
Kpl = 5 
Kp2 = 5 
Kx5 = 10 

# No feed-back 
alphal = 0 
alpha2 = 0 

# Negative feed-back from X4 to 

# alphal = 50 

# alpha2 = 0 

# InitVar 
xl = 1.0 
x2 = 1.0 
x3 = 1.0 
x4 = 1.0 



G Gene Regulatory Networks 
G.l Two-enzyme GRN 

FIX: Wl Gl G2 

# Reactions 
Rl: 

Gl > El 

Vmaxl*(l/(l+kl*El**hl))*(l/(l+k2*P**li2))*(l/(l+k3*S**h3)) 

R2: 

G2 > E2 

Vmax2*(l/(l+k4*E2**h4))*(l/(l+k5*El**h5)) 

R3: 

El > Wl 
E1*0. 650524 

R4: 

E2 > Wl 
E2*0. 223366 

# Fixed Species 
Wl = 1 

Gl = 1 
G2 = 1 

# Variable Species 
El = 1 

E2 =1 
S = 1 
P = 1 

# Parameters 
Vmaxl = 3.75422 
kl = 6.48859 

hi = 4.30682 
k2 = 0.402033 
h2 = 1.00019 
k3 = 0.43719 
h3 = 0.869487 
Vmax2 = 2.46437 
k4 = 0.86132 
h4 = 3.6914 
k5 = 4.90228 
h.5 = 3.14973 

G.2 Three-enzyme GRN 

FIX: Wl Gl G2 G3 

# Reactions 
Rl: 

Gl > El 

Vmaxl* ((l+kl*P**hl)*(l+k3*El**h3)-l)/((l+kl*P**hl)*(l+k2*S**h2)*(l+k3*El**h3)) 

R2: 

G2 > E2 
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Vmax2* ( (l+k4*S**h4) -1) / ( (l+k4*S**h4) * (l+k5*P**h5) ) 

R3: 

G3 > E3 

Vmax3* ( (l+k6*S**h6) * (l+k7*P**h7) * (l+k8*E3**h8) -1) / 
( (l+k6*S**h6) * (l+k7*P**h7) * (l+k8*E3**h8) ) 

R4: 

El > Wl 
El*l. 09882 

R5: 

E2 > Wl 
E2*l. 27824 

R6: 

E3 > Wl 
E3*l. 29683 

# Fixed Species 
Wl = 1 

Gl = 1 
G2 = 1 
G3 = 1 

# Variable Species 
El = 1 

E2 = 1 
E3 = 1 
S = 1 
P = 1 

# Parameters 
Vmaxl = 2.83741 
kl = 1.63564 

hi = 1.49258 
k2 = 2.00275 
b.2 = 0.247268 
k3 = 8.92731 
h3 = 0.157243 

Vmax2 = 3.28334 
k4 = 1.01241 
h.4 = 0.187473 
k5 = 0.0641848 
b.5 = 0.487284 

Vmax3 = 1.24576 
k6 = 8.5334 
h.6 = 0.216253 
k7 = 9.47921 
h.7 = 1.63791 
k8 = 2.45668 
h.8 = 3.99067 
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G.3 Five-enzyme GRN 

FIX: Wl Gl G2 G3 G4 G5 

# Reactions 
Rl: 

Gl > El 

Vmaxl* ( (l+k2*El**h2) -1) / ( (l+kl*E5**hl) * (l+k2*El**h2) ) 

R2: 

G2 > E2 

Vmax2* ( (l+k3*P**h3) * (l+k4*S**h4) -1) / ( (l+k3*P**h3) * (l+k4*S**h4) ) 

R3: 

G3 > E3 

Vmax3* ( (l+k5*E5**h5) -1) / ( (l+k5*E5**h5) * (l+k6*E3**h6) ) 

R4: 

G4 > E4 

Vmax4* ( (l+k7*E5**h7) -1) / ( (l+k7*E5**h7) ) 

R5: 

G5 > E5 

Vmax5*((l+k8*S**h8)*(l+k9*P**h9)-l)/((l+k8*S**8)*(l+k9*P**h9)*(l+kl0*E5**hl0)) 

R6: 

El > Wl 
E1*0. 516899 

R7: 

E2 > Wl 
E2*0. 867095 

R8: 

E3 > Wl 
E3*l. 41811 

R9: 

E4 > Wl 
E4*0. 797611 

RIO: 

E5 > Wl 
E5*0. 256305 

# Fixed Species 
Wl = 1 

Gl = 1 
G2 = 1 
G3 = 1 
G4 = 1 
G5 = 1 

# Variable Species 
El = 1 

E2 = 1 
E3 = 1 
E4 = 1 
E5 = 1 
S = 1 
P = 1 
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# Parameters 
Vmaxl = 2.43445 
kl = 6.6671 
hi = 1.54412 
k2 = 9.67675 
h2 = 2.53249 
Vmax2 = 1.24985 
k3 = 2.10143 
h3 = 0.0121673 
k4 = 1.97132 
h4 = 1.04927 
Vmax3 = 2.30057 
k5 = 6.64561 
h5 = 1.71241 
k6 = 0.433479 
h6 = 2.48621 
Vmax4 = 1.34903 
k7 = 3.16177 
h7 = 1.59593 
Vmax5 = 0.183229 
k8 = 0.896061 
h8 = 0.580101 
k9 = 0.0406431 
h9 = 0.612953 
klO = 0.085952 
hlO = 0.938131 



G.4 Six-enzyme GRN 

FIX: Wl Gl G2 G3 G4 G5 G6 S P 

# Reactions 
Rl: 

Gl > El 

Vmaxl*((l+kl*El**hl)-l)/((l+kl*El**hl)*(l+k2*E6**h2)) 

R2: 

G2 > E2 

Vmax2* ( (l+k3*E4**h3) -1) / ( (l+k3*E4**h3) ) 

R3: 

Gl > E3 

Vmax3* (1/ (l+k5*P**h5) ) * ( (k6*E6**h6) / (l+k6*E6**h6) ) 

R4: 

Gl > E4 

Vmax4* ( (k7*E6**h7) / (l+k7*E6**h7) ) 

R5: 

Gl > E5 

Vmax5* ( (k9*E6**h9) / (l+k9*E6**h9) ) 

R6: 

Gl > E6 

Vmax6*((l+kll*S**hll)*(l+kl2*P**hl2)-l)/((l+kll*S**lill)*(l+kl2*P**hl2)) 

R7: 

El > Wl 
E1*0. 186924 

R8: 

E2 > Wl 
E2*l. 88726 

R9: 

E3 > Wl 
E3*l. 57073 

RIO: 

E4 > Wl 
E4*1.40082 

Rll: 

E5 > Wl 
E5*l. 67811 

R12: 

E6 > Wl 
E6*0. 440807 

# Fixed Species 
S = 100 

P = 1 
Wl = 1 
Gl = 1 
G2 = 1 
G3 = 1 
G4 = 1 
G5 = 1 
G6 = 1 
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# Variable Species 

El = 1 

E2 = 1 

E3 = 1 

E4 = 1 

E5 = 1 

E6 = 1 



# Parameters 
Vmaxl = 1.4411 
kl = 2.40536 
hi = 0.152005 
k2 = 8.1128 
b.2 = 1.40023 
Vmax2 = 3.06359 
k3 = 7.698 
h3 = 1.65138 
Vmax3 = 2.01559 
k5 = 0.0272647 
h.5 = 0.106304 
k6 = 5.03767 
h.6 = 0.994795 
Vmax4 = 2.07882 
k7 = 8.22248 
b.7 = 1.56628 
Vmax5 = 2.15233 
k9 = 3.33499 
h.9 = 1.44571 
Vmax6 = 0.244187 
kll = 0.77059 
hll = 0.635955 
kl2 = 0.028323 
hl2 = 1.19337 



